#ifndef PROTEIN_GRAPH_PARSING_H_
#define PROTEIN_GRAPH_PARSING_H_

#include "Graph.h"

namespace masugo
{

class protein_graph_parsing
{
public:
	typedef Graph<std::string, std::string, int> FGraph;

	static std::deque<FGraph> read_graphs(const char* filename)
	{
		try
		{
			std::ifstream input(filename, std::ios::in);
			std::deque<FGraph> graphs;
			while(input.good())
			{
				std::string graph_name;
				std::getline(input, graph_name);
				if(!graph_name.empty())
				{
					FGraph graph=read_graph(input);
					graph.set_name(graph_name);
					graphs.push_back(graph);
				}
			}
			return graphs;
		}
		catch(const std::exception& e)
		{
			throw std::runtime_error(utilities::ParametricString("Failed to read protein graphs, invalid format for the file '{}'").arg(filename).str());
		}
	}

private:
	static std::string residue_id(const utilities::Table& atoms_table, const int row)
	{
		return (atoms_table.value<std::string>(row, "Chain-name")+atoms_table.value<std::string>(row, "Residue-serial"));
	}

	static bool is_side_chain(const utilities::Table& atoms_table, const int row)
	{
		const std::string atom_name=atoms_table.value<std::string>(row, "Atom-name");
		return (atom_name!="C" && atom_name!="CA" && atom_name!="N" && atom_name!="O" && atom_name!="OXT");
	}

	static bool are_chain_neighbours(const utilities::Table& atoms_table, const int row_a, const int row_b)
	{
		return (atoms_table.value<std::string>(row_a, "Chain-name")==atoms_table.value<std::string>(row_b, "Chain-name")
				&& abs(atoms_table.value<int>(row_a, "Residue-serial")-atoms_table.value<int>(row_b, "Residue-serial"))==1);
	}

	static FGraph read_graph(std::istream& input)
	{
		FGraph graph;

		const utilities::Table atoms_table(input);
		const utilities::Table contacts_table(input);

		for(int i=0;i<atoms_table.rows_count();i++)
		{
			if(atoms_table.value<int>(i, "ID")!=i) { throw std::invalid_argument(utilities::ParametricString("Function {} encountered invalid atoms table").arg(__FUNCTION__).str()); }
			graph.set_vertex_value(residue_id(atoms_table, i), atoms_table.value<std::string>(i, "Residue-name"));
		}

		for(int i=0;i<contacts_table.rows_count();i++)
		{
			const int a_i=contacts_table.value<int>(i, "A");
			const int b_i=contacts_table.value<int>(i, "B");
			if(a_i!=b_i && is_side_chain(atoms_table, a_i) && is_side_chain(atoms_table, b_i))
			{
				const std::string a_residue_id=residue_id(atoms_table, a_i);
				const std::string b_residue_id=residue_id(atoms_table, b_i);
				if(a_residue_id!=b_residue_id)
				{
					graph.connect(a_residue_id, b_residue_id, are_chain_neighbours(atoms_table, a_i, b_i) ? 0 : 1);
				}
			}
		}

		return graph;
	}
};

}

#endif /* PROTEIN_GRAPH_PARSING_H_ */
